{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###### Content under Creative Commons Attribution license CC-BY 4.0, code under BSD 3-Clause License © 2018  by D. Koehn, notebook style sheet by L.A. Barba, N.C. Clementi"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<link href=\"https://fonts.googleapis.com/css?family=Merriweather:300,300i,400,400i,700,700i,900,900i\" rel='stylesheet' >\n",
       "<link href=\"https://fonts.googleapis.com/css?family=Source+Sans+Pro:300,300i,400,400i,700,700i\" rel='stylesheet' >\n",
       "<link href='http://fonts.googleapis.com/css?family=Source+Code+Pro:300,400' rel='stylesheet' >\n",
       "<style>\n",
       "\n",
       "@font-face {\n",
       "    font-family: \"Computer Modern\";\n",
       "    src: url('http://mirrors.ctan.org/fonts/cm-unicode/fonts/otf/cmunss.otf');\n",
       "}\n",
       "\n",
       "\n",
       "#notebook_panel { /* main background */\n",
       "    background: rgb(245,245,245);\n",
       "}\n",
       "\n",
       "div.cell { /* set cell width */\n",
       "    width: 800px;\n",
       "}\n",
       "\n",
       "div #notebook { /* centre the content */\n",
       "    background: #fff; /* white background for content */\n",
       "    width: 1000px;\n",
       "    margin: auto;\n",
       "    padding-left: 0em;\n",
       "}\n",
       "\n",
       "#notebook li { /* More space between bullet points */\n",
       "margin-top:0.5em;\n",
       "}\n",
       "\n",
       "/* draw border around running cells */\n",
       "div.cell.border-box-sizing.code_cell.running { \n",
       "    border: 1px solid #111;\n",
       "}\n",
       "\n",
       "/* Put a solid color box around each cell and its output, visually linking them*/\n",
       "div.cell.code_cell {\n",
       "    background-color: rgb(256,256,256); \n",
       "    border-radius: 0px; \n",
       "    padding: 0.5em;\n",
       "    margin-left:1em;\n",
       "    margin-top: 1em;\n",
       "}\n",
       "\n",
       "\n",
       "div.text_cell_render{\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    line-height: 140%;\n",
       "    font-size: 110%;\n",
       "    width:680px;\n",
       "    margin-left:auto;\n",
       "    margin-right:auto;\n",
       "}\n",
       "\n",
       "/* Formatting for header cells */\n",
       ".text_cell_render h1 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-style:regular;\n",
       "    font-weight: bold;    \n",
       "    font-size: 250%;\n",
       "    line-height: 100%;\n",
       "    color: #004065;\n",
       "    margin-bottom: 1em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       ".text_cell_render h2 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: bold; \n",
       "    font-size: 180%;\n",
       "    line-height: 100%;\n",
       "    color: #0096d6;\n",
       "    margin-bottom: 0.5em;\n",
       "    margin-top: 0.5em;\n",
       "    display: block;\n",
       "}\t\n",
       "\n",
       ".text_cell_render h3 {\n",
       "    font-family: 'Merriweather', serif;\n",
       "\tfont-size: 150%;\n",
       "    margin-top:12px;\n",
       "    margin-bottom: 3px;\n",
       "    font-style: regular;\n",
       "    color: #008367;\n",
       "}\n",
       "\n",
       ".text_cell_render h4 {    /*Use this for captions*/\n",
       "    font-family: 'Merriweather', serif;\n",
       "    font-weight: 300; \n",
       "    font-size: 100%;\n",
       "    line-height: 120%;\n",
       "    text-align: left;\n",
       "    width:500px;\n",
       "    margin-top: 1em;\n",
       "    margin-bottom: 2em;\n",
       "    margin-left: 80pt;\n",
       "    font-style: regular;\n",
       "}\n",
       "\n",
       ".text_cell_render h5 {  /*Use this for small titles*/\n",
       "    font-family: 'Source Sans Pro', sans-serif;\n",
       "    font-weight: regular;\n",
       "    font-size: 130%;\n",
       "    color: #e31937;\n",
       "    font-style: italic;\n",
       "    margin-bottom: .5em;\n",
       "    margin-top: 1em;\n",
       "    display: block;\n",
       "}\n",
       "\n",
       ".text_cell_render h6 { /*use this for copyright note*/\n",
       "    font-family: 'Source Code Pro', sans-serif;\n",
       "    font-weight: 300;\n",
       "    font-size: 9pt;\n",
       "    line-height: 100%;\n",
       "    color: grey;\n",
       "    margin-bottom: 1px;\n",
       "    margin-top: 1px;\n",
       "}\n",
       "\n",
       "    .CodeMirror{\n",
       "            font-family: \"Source Code Pro\";\n",
       "\t\t\tfont-size: 90%;\n",
       "    }\n",
       "/*    .prompt{\n",
       "        display: None;\n",
       "    }*/\n",
       "\t\n",
       "    \n",
       "    .warning{\n",
       "        color: rgb( 240, 20, 20 )\n",
       "        }  \n",
       "</style>\n",
       "<script>\n",
       "    MathJax.Hub.Config({\n",
       "                        TeX: {\n",
       "                           extensions: [\"AMSmath.js\"], \n",
       "                           equationNumbers: { autoNumber: \"AMS\", useLabelIds: true}\n",
       "                           },\n",
       "                tex2jax: {\n",
       "                    inlineMath: [ ['$','$'], [\"\\\\(\",\"\\\\)\"] ],\n",
       "                    displayMath: [ ['$$','$$'], [\"\\\\[\",\"\\\\]\"] ]\n",
       "                },\n",
       "                displayAlign: 'center', // Change this to 'center' to center equations.\n",
       "                \"HTML-CSS\": {\n",
       "                    styles: {'.MathJax_Display': {\"margin\": 4}}\n",
       "                }\n",
       "        });\n",
       "</script>\n"
      ],
      "text/plain": [
       "<IPython.core.display.HTML object>"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Execute this cell to load the notebook's style sheet, then ignore it\n",
    "from IPython.core.display import HTML\n",
    "css_file = '../style/custom.css'\n",
    "HTML(open(css_file, \"r\").read())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Accuracy of Finite-difference operators\n",
    "\n",
    "In the [FD introduction lesson](http://nbviewer.jupyter.org/github/daniel-koehn/Theory-of-seismic-waves-II/blob/master/03_Intro_finite_differences/1_fd_intro.ipynb), we introduced finite-difference approximations of the first derivative by neglecting the limit $dx \\rightarrow 0$, which is an approximation of the tangent slope of $f(x)$ at the point $x$ by the slope of the secant through two points in the neighborhood. In order to characterize the accuracy of FD operators, I introduce a more systematic approach by Taylor series expansion."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Taylor series and accuracy of FD operators\n",
    "\n",
    "We derive the **forward FD operator**:\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^+ \\approx \\frac{f(x+dx)-f(x)}{dx}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "by calculating the **Taylor series expansion** of $f(x+dx)$ around the point $x$:\n",
    "\n",
    "\\begin{equation}\n",
    "f(x+dx) \\approx f(x) + \\frac{\\partial f(x)}{\\partial x} dx + \\frac{1}{2}\\frac{\\partial^2 f(x)}{\\partial x^2} dx^2\\; +\\; ...\\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "After rearranging to $\\frac{\\partial f(x)}{\\partial x}$, we get\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^+ \\approx \\frac{f(x+dx)-f(x)}{dx} {\\color{red}{-\\frac{1}{2}\\frac{\\partial^2 f(x)}{\\partial x^2} dx}}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "Because the remainder term of the Taylor series (marked in red), has a term of order dx,  the FD operator\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^+ \\approx \\frac{f(x+dx)-f(x)}{dx}+{\\color{red}{\\mathcal{O}(dx)}} \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "is a **first order FD operator**. The term first order refers to the order of the Taylor series remainder and should not be confused by the order of a derivation. The error of a first order FD operator decreases linear with $dx$.\n",
    "\n",
    "Using the Taylor series expansion of \n",
    "\n",
    "\\begin{equation}\n",
    "f(x-dx) \\approx f(x) - \\frac{\\partial f(x)}{\\partial x} dx + \\frac{1}{2}\\frac{\\partial^2 f(x)}{\\partial x^2} dx^2\\; -\\; ...\\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "we can show that the **backward FD operator**\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^- \\approx \\frac{f(x)-f(x-dx)}{dx}+{\\color{red}{\\mathcal{O}(dx)}} \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "is also a first order operator. \n",
    "\n",
    "What is the accuracy of the **central FD operator**?\n",
    "To answer this question, we expand the Taylor series of $f(x+dx)$ and $f(x-dx)$ by an additional term:\n",
    "\n",
    "\\begin{align}\n",
    "f(x+dx) &\\approx f(x) + \\frac{\\partial f(x)}{\\partial x} dx + \\frac{1}{2}\\frac{\\partial^2 f(x)}{\\partial x^2} dx^2 + \\frac{1}{6}\\frac{\\partial^3 f(x)}{\\partial x^3} dx^3 + ...\\nonumber\\\\\n",
    "f(x-dx) &\\approx f(x) - \\frac{\\partial f(x)}{\\partial x} dx + \\frac{1}{2}\\frac{\\partial^2 f(x)}{\\partial x^2} dx^2 - \\frac{1}{6}\\frac{\\partial^3 f(x)}{\\partial x^3} dx^3 + ...\\nonumber\n",
    "\\end{align}\n",
    "\n",
    "Subtracting $f(x-dx)$ from $f(x+dx)$ leads to \n",
    "\n",
    "\\begin{equation}\n",
    "f(x+dx) - f(x-dx) \\approx 2\\frac{\\partial f(x)}{\\partial x} dx + \\frac{2}{6}\\frac{\\partial^3 f(x)}{\\partial x^3} dx^3 + ...\\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "After solving for $\\frac{\\partial f(x)}{\\partial x}$, we derived the central FD operator\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^c \\approx \\frac{f(x+dx)-f(x-dx)}{2dx} {\\color{red}{- \\frac{1}{6}\\frac{\\partial^3 f(x)}{\\partial x^3} dx^2}}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "In contrast to the forward and backward operator, the Taylor series remainder of the central operator (marked in red), has an order $dx^2$\n",
    "\n",
    "\\begin{equation}\n",
    "\\biggl(\\frac{\\partial f(x)}{\\partial x}\\biggr)^c \\approx \\frac{f(x+dx)-f(x-dx)}{2dx} + {\\color{red}{\\mathcal{O}(dx^2)}}. \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "Therefore, the central approximation is a **second order FD operator**, meaning that the error decreases quadratically with $dx^2$. This explains why the central FD operator is more accurate than the forward and backward operators, as we have seen in the FD introductory lesson, where we compared the first derivation of the Gaussian using the different FD operators with the analytical solution.\n",
    "\n",
    "<img src=\"images/FD_error_1st_derivative_gauss.png\" width=\"95%\">\n",
    "\n",
    "##### Bonus Exercise\n",
    "\n",
    "Quantify the errors for the first derivative to the Gaussian approximated by the forward, backward and central operators using the remainder of the Taylor series expansion. Extend the notebook [FD Introduction](https://nbviewer.org/github/daniel-koehn/Theory-of-seismic-waves-II/blob/master/03_Intro_finite_differences/1_fd_intro.ipynb) to plot the remainder of the Taylor series expansion in comarison with the errors between the numerical and analytical solutions.\n",
    "\n",
    "##### Exercise\n",
    "\n",
    "What is the accuracy of the 3-point FD operator\n",
    "\n",
    "\\begin{equation}\n",
    "\\frac{\\partial^2 f(x)}{\\partial x^2} \\approx \\frac{f(x+dx)-2f(x)+f(x-dx)}{dx^2} \\nonumber\n",
    "\\end{equation}\n",
    "\n",
    "which approximates the 2nd derivative of $f(x)$?\n",
    "\n",
    "Hint: Use the derivation of the 3-point FD operator from the [FD introduction lesson](http://nbviewer.jupyter.org/github/daniel-koehn/Theory-of-seismic-waves-II/blob/master/03_Intro_finite_differences/1_fd_intro.ipynb)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## We learned:\n",
    "\n",
    "* How to estimate the accuracy of FD operators by Taylor series expansion"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
